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COMPUTER PROGRAM FOR EVALUATION OF BLQCH-GRUENEISEN PARAMETERS 
OF METALS AND EVALUATION OF ELECTRICAL RESISTIVITY 
OF TANTALUM AS A FUNCTION OF TEMPERATURE 
by Thor T. Semler and John P. Riehl 
Lewis Research Center 

SUMMARY 


A computer program using nonlinear functional minimization has been written to ob- 
tain least squares solutions, from experimental measurements, for the constant terms 
of the Bloch-Grueneisen relation 


Prjl(T) - A 


/> e R /T 

/ T\ 5 j x 5 dx 

V R/ J Q (e x - 1)(1 - e" x ) 


where p T (T) is the "ideal" electrical resistivity of metals as a function of temperature, 
A is a constant of the metal, T is the temperature in K, and 0^ is the characteristic 
temperature of resistivity. 

The metal tantalum has been analyzed by using the code, and a typical result is 
= 217. 54 K and A = 39. 95 microohm-centimeters between 10 and 250 K. 


INTRODUCTION 

A fundamental constant in Ohm's law is the electric resistance of the conductor R. 

For a particular conductor of length l (in cm) and uniform cross-sectional area a (in 
2 

cm ), R (in ohms) may be computed as 


R = 


pl 


a 


(1) 



where p is the electrical resistivity in ohm -centimeters. As the electrical resistance 
is a quantity of great interest in both engineering and solid-state physics (refs. 1 to 3), 
it is important to be able to determine the electrical resistivity and then compute the re- 
sistance of a conductor. 

The resistivity of a metal is a function of temperature. On approaching very low 
temperatures, near absolute zero, the electrical resistivity assumes a constant value 
(neglecting the region in which some metals become superconductors) p , called the re- 
sidual resistivity. This quantity p Q arises from imperfections, impurities, and strains 
in the metal lattice and must be determined for each individual sample. 

The total resistivity p may be divided into two portions, the residual resistivity and 
the temperature -dependent resistivity p^( T): 

P = P 0 + P T (T) (2) 

This division is known as Matthies sen's rule (refs. 3 and 4). 

It is possible to derive a formula for the temperature-dependent resistivity p T (T) 
over a large temperature range from certain approximations about the interactions of 
conduction electrons and the metallic lattice vibrations (refs. 5 to 8). The formula (3) 
so derived is referred to as the Bloch-Grueneisen relation 

/ flp/T 

(3) 

(e x - 1)(1 - e _x ) 

where A is a constant of the metal and 0^ is the characteristic temperature of resis- 
tivity. The Bloch-Grueneisen relation is widely applied because it provides a good ap- 
proximation for the temperature-dependent resistivity for many metals. 

Because of the difficulties in the form of the relation, many rule-of-thumb techniques 
have been evolved by experimenters to evaluate A and 0 R from experimental data. Un- 
fortunately, many of these rule-of-thumb techniques are rather crude (ref. 4). 

The computer program described in the section ANALYSIS allows one to compute in 
the least -squares sense the best values of A and from all the experimental data one 
might wish to use. This computer program is then used to obtain values of A and 0 R 
for metallic tantalum. 



ANALYSIS 


Given values of p T that have been measured as a function of temperature, one 

would like the values of 0 R and A which are ’best" in the least-squares sense. This 

means that one must form the function f(A, 9 ) - Y, , , , 71 2 

\ y rT,i, measured ^T, 1, calculated 

and minimize it by the variation of 0 R and A. - 1 

By the rules of ordinary calculus, f(A, 9) without constraint obtains a local minimum 
or reaches a saddle point when the gradient Vf(A, 9 ) = 0 at particular values of A and 
9. For those functions Z(x,p) that are linear in p the gradient requirement produces a 
set of simultaneous linear equations. But for functions that are nonlinear in p, this re- 
quirement is not easily met. One is confronted with a set of nonlinear simultaneous 
equations in p to be solved. Such is the case for the Bloch-Grueneisen equation. Its 
gradient has components 




where NDP is the number of data points. 
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The second component was arrived at by using the chain rule and Leibnitz's rule, which 


if F(x) = ^ (x) 


f(x,y)dy is continuously differentiable, then 


F 



— f(x,y)dy - f(x , cp(x)) 
0X 


cp'(x) + f(x,if/(x)) ■ <p'(x) 


To accomplish the unconstrained function minimization, the method of Fletcher and 
Powell was employed (ref. 9). The Fletcher and Powell algorithm is a modification of a 
method of Davidon (ref. 10). It is a powerful and general method for finding the local 
minimum of a general function f(x). 

Central to the method is a symmetric positive definite matrix jf- which is updated 
at each iteration i. The current direction of motion is supplied by when it is 
multiplied with the current change gradient vector. An iteration is described by the 
following: If is any positive definite matrix, usually the identity matrix /, on the 

first iteration only, then 


§1 = 

Choose a = a. by minimizing f(x^ + aST^); this straight line minimization is done 
with cubic interpolation: 


= “ s i 


x. , = x. + a. 
l+l i i 


■*i+l - 


where the matrices 


si ^ and are defined by 




and 


y. being the transpose of y.. 


a. 





The numerators of d i and #. are both matrices, while the denominators are scalars. 
Fletcher and Powell (ref. 9) prove the following: 

(1) The matrix is positive and definite for all i. As a consequence, the method 
will always converge since 

— f(x. + aS.) | Q = -Vf(x-) Jfj Vf(x.)< 0 
d(v 

That is, the function f is initially decreasing along the direction S.. So that the func- 
tion can be decreased at each iteration by minimizing along S- . 

(2) When the method is applied to the quadratic matrix equation q(x) = a + bx + xVx 
and x is a vector of n dimensions , 

(a) The directions S. (or equivalently a.) are <4 conjugate, that is, 

S^ Sj = 0 for i 4 j . This condition leads to a minimum in n steps . 

(b) The matrix ^ converges to the inverse of the Hessian, that is, the matrix 

of second partial derivatives after n iterations, When applied to a gen- 

eral function f(x), tends to the inverse of the Hessian evaluated at the minimum. 
The Fletcher -Powell algorithm is represented by the flow chart in figure 1. 


CALCULATION OF BLOCH-GRUENEISEN RELATION 

The integral portion of the Bloch-Grueneisen relation was calculated by using a mod- 
ified Simpson’s rule integration scheme. This scheme, programmed as subprogram 
SIMPS 1, adapts to regions where more points are required to obtain an accurate result. 
Had the integral been too expensive (in terms of computer time) to compute, a spline ap- 
proximation to tabular results of the integral could have been made. It was our experi- 
ence that the Bloch-Grueneisen relation and the ensuing least -squares function could be 
calculated in very little time by using the modified Simpson's rule routine. 

In the program the exponent 5 of equation (3) may be varied, as indicated in the com- 
ments card. This permits the user to use other than a fifth-order Bloch-Grueneisen 
relation. Input is described in appendix A. A flow chart of the main program is shown 
in figure 2. A listing of the program is given in appendix B. 


DATA USED IN EVALUATION OF TANTALUM 


Experimental measurements of the electrical resistivity of the metal tantalum have 
been analyzed by using this code. Only experimental values of the electrical resistivity 
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in the temperature region from above 0 K (actually 10 K as tantalum is a superconductor) 
to about 400 K have been used in the program. The values have been taken from this re- 
gion since the least -squares fit of the parameters is relatively insensitive to data outside 
the range 0 K to a few times 0-^ K. 

Cox Data in Temperature Range 77 to 373 K 

Cox in 1943 performed a series of experiments to determine both the thermal and 
electrical conductivities of tungsten and tantalum (ref. 11). This series of experiments 
yielded three values of the electrical resistivity of tantalum. 

The sample used was a tantalum wire about 40 centimeters long and about 0. 0254 
centimeter in diameter. The wire was aged by passing as high a current as possible 
through it without evaporating it. The tantalum sample was aged at both 1800° and 
2000° C for a total time of 2750 hours. The resistance at zero power input was plotted 
at a function of aging time. The resistance decreased rapidly at first and finally reached 
a constant value; at that point, aging was ceased. The chemical purity quoted for the 
sample was 99.9 percent. 

After aging, the sample was immersed in baths of boiling liquid nitrogen, ice water, 
and boiling water, readings of voltage and current across the sample were taken, and the 
resistivity was computed; these data are given in table I. 


White and Woods Data in Temperature Range 10 to 295 K 

White and Woods, in a series of experiments to determine the electrical and thermal 
resistivity of the transition elements, report 21 values of the ideal resistivity of tantalum 
(ref. 12). These results were obtained by subtracting the residual resistivity from the 
total resistivity of the sample at a temperature, and they are shown in table II. 

The specimen was mounted in a cryostat. One end of the specimen was soldered to a 
copper post, and a specimen heater was attached to the other end. Copper wires were 
attached to intermediate points of the rod to act as electrical potential leads for the resis- 
tivity measurements. The specimen was surrounded by low-pressure helium gas to pre- 
serve temperature equilibrium in the cryostat. 

The purities of the three samples of tantalum used are quoted as 99. 9 percent, 

"high", and 99.9 percent. All samples were vacuum annealed to remove as much work 
hardening as possible. 

White and Woods suggest that the electrical resistivity of tantalum at lower temper - 

3 8 5 

atures follows more nearly a T ' proportionality than a T proportionality. However, 



in the section of their report devoted to the error analysis of the electrical resistivity 
measurements they indicate the difficulties in determining the low-temperature ideal 
resistivities . 


RESULTS AND DISCUSSION 


The program was executed by using 19 of White and Woods experimental values. 

This involved a temperature range from 10 to 250 K. The resultant values of the param- 
eters were A = 39. 95 microohm -centimeters and 0^ = 217. 54 K. The resultant fit of 
the calculated data to the experimental points is shown in figure 3. It can be seen that 
the agreement is excellent. The tabular results are shown in appendix C. 

The three values of Cox were analyzed by using the code. This involved a temper- 
ature range of 77. 3 to 373.4 K. The result of this analysis were A = 39. 51 microohm- 
centimeters and = 210. 77 K. The resulting fit of the data is in figure 4 . 

Because of the rather limited nature of the Cox data, both the White and Woods data 
and the Cox data were analyzed together. The resulting fit is shown in figure 5. It can 
be clearly seen that the Cox data appear higher than the White and Woods data. This in- 
dicates either that the residual resistivity of the Cox sample had not been subtracted 
from the individual values or that the sample had been insufficiently annealed. Thus, the 
Cox data have been rejected in following the analysis. 

As an illustrative example of the utility of the code, the White and Woods data have 
been analyzed parametrically by using the highest temperature involved as a parameter. 
The lowest temperature in all these cases was 10 K. The results are shown in figures 8 
and 7. The results are tabulated in table HI. They show both A and 0^ increasing to 
their asymptotic values. These are typical results for metallic samples (ref. 4). The 
characteristic temperature of resistivity 0-^ is not to be confused with 0 D , the Debye 
characteristic temperature (ref. 13). 

It should be indicated at this point that while these results might be obtained by ex- 
tensive hand calculation, the results shown were obtained in a fraction of a minute by an 
IBM 7094 -13 computer. 


SUMMARY OF RESULTS 

A computer program was written for evaluation of Bloch-Grueneisen parameters of 
metals and evaluation of electrical resistivity of tantalum as a function of temperature, 
and the following results were obtained: 



1. The values of the Bloch-Grueneisen parameters for the data of White and Woods 
were a constant of the metal A = 39. 95 microohm-centimeters and a characteristic 
temperature of resistivity Q ^ = 217. 54 K over the temperature range from 10 to 250 K. 

2. The results of Cox apparently were not adjusted to ideal resistivity values. 

3. When the Bloch-Grueneisen parameters were plotted as a function of the highest 
temperature used, for tantalum, they reached their asymptotic values at roughly a tem- 
perature of 

Lewis Research Center, 

National Aeronautics and Space Administration, 

Cleveland, Ohio, April 19, 1971, 

129-02. 



APPENDIX A 


PROGRAM DESCRIPTION 
Program Input Data 

The card input to the program consists of a temperature T., ^ , the source of the 

data, and the first guess for A and The format for these cards’ is 2F10. 0,3X, A6. 
The user supplies as many cards with T. , p T and the source as needed. The last 
card contains the first guess for A and 0 ^ with the field for the source left blank. 

This card causes the start of the least -squares curve fit. If the user wants the proceed- 
ing values of A and 0 -^ as first guesses, the last card should be left entirely blank. 

The user may add coding of his own, starting after statement number 7 (card num- 
ber 37). By punching END on a data card starting in card column 24, the user may exe- 
cute this coding. This END card causes the program to transfer to statement 7. 


Program Output 

The program output consists of a tabulation of the input data along with p T, , the 
calculated value of p T ^ as defined by the Bloch-Grueneisen relation, and the difference 
between p c Tj and p T .. The values of A, 0 R , and the sum of the differences squared 
are also printed. An example of this output is found in appendix C. 



APPENDIX B 


FORTRAN LISTING OF PROGRAM AND OUTPUT 


COMMON /BLOCK/ T ( luG ) , RHC 1 100 ) * RHCCt 100 ) , I A I 

COMMON /SPACE/ WORK (10 ) A 2 

COMMON /ESTIM/ EST.EPS.L1 MIT ,IER A 3 

EXTERNAL FAT A 4 

DIMENSION GRA0I2I , XI2I , SCURCE(IOO) A 5 

EQUIVALENCE (XU), A), (X<2). THETA) A 6 

DATA BLANK/1H / A 7 

DATA ENU/3HE ND/ A 8 

DATA N/2/ A 9 

WRITE 16,8) A 10 

1 I=i A 11 

2 READ (3,12) T (I ) , RH C ( I ) , SOURCE ( I ) A 12 

IF ( SOURCE (I J-BLANK) 3,4,3 A 13 

3 CONTINUE A 14 

If ( SOURCE (I ) .EQ. END) GO TO 7 A 15 

1=1 + 1 A 16 

GO TO 2 A 17 

4 LON TI NCE A 1 8 

IF (TU).EQ.O.) GO TO 5 A 19 

C MAKE A FIRST GUESS AT A AND THETA A 20 

A =T ( I ) A 21 

THE TA =RHO (I ) A 22 

3 CONTINUE A 23 

1=1-1 A 24 

C CALL THE FLETCHER - POWELL SUBROUTINE A 25 

CALL FLTPWL ( FA T , N , X , VA L , GRAD ) A 26 

IF ( IER.NE.O) WRITE (6,9) IER A 27 

WRITE (6,11) A 28 

00 6 J=1,I A 29 

UIF=RHO( JI-RHOC ( J) A 30 

WRITE (6,13) SOURCE ( J) , T( J) ,RHO( J) ,RHOC(J ) ,DIF A 31 

6 CONTINUE A 32 

WRITE (6,10) A 33 

wRITt (6,14) A , THE TA , VAL A 34 

WRITE (6,8) A 35 

GO TO 1 A 36 

7 CONTINUE A 37 

C PEKFURM ANY OTHER CALCULATIONS HERE A 38 

STOP A 39 

C A 40 

C A 41 

8 FORMAT ( iHl ) A 42 

9 FORMAT (10X,4HIER=,I2) A 43 

10 FORMAT { 1HK,24X ,1HA ,16 X ,5HTHET A,8X,23HSUM OF ( DIFFERENCES )**2) A 44 

11 FORMAT 11HK,16X,6HS0URCE,17X,11HTEMPERATURE»13X»3HRH0,12X» 14HRH0 C A 45 

1ALCULATED , 8X ,1 OHDI FF E RE NCE) A 46 

12 FORMAT ( 2Flu« 0 , 3 X , A6 ) A 47 

13 FORMAT (17X,A6,7X,4F20.6) A 48 

14 FORMAT (10X.3F20.6) A 49 

END A 50- 
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J o u 


SUBROUTINE FAT I N , Q , VLL ,GRDD ) B 1 

THIS SUBROUTINE CALCULATES THE LEAST SQUARES FUNCTION FIX I B 2 

ANU THE GRAD I ENT OF THE SAME B 3 

IN THIS PART XII) =A , X(2)=THETA B A 

COMMON /BLOCK/ T I IOC ) , RHC ( 100 ) , BGRI 100 ) , NCAS ES B 5 

COMMON /EXPK/ K B 6 

DIMENSION X(2), GRDD(l), Q(l), GRADI2) B 7 

EXTERNAL FUNKY 6 8 

DO 1 1=1, N B 9 

1 XII) =Q (II B 1C 

JOKE =0 B 11 

GO TO 2 B 12 

ENTRY BLAST! Z,Y, VLL) B 13 

X(1)=Z B 14 

X ( 2 ) = Y B 15 

JOKE =1 B 16 

2 CONTINUE B IT 

VAL =0. B 18 

DO 3 1=1,2 B 19 

3 GRAD (I ) =0. B 20 

AK=K B 37 

DO 4 1=1, NCA SE S B 21 

C CALCULATE THE BLOCH - GRUENEISEN RELATIONSHIP B 22 

X2DTI=X(2)/T(I) B 23 

TIDX2=(T(I)/X(2) )**K B 24 

XX = S1MPS1 (0. ,X2DTI , FUNKY, U B 25 

BGRII )=X(1I*TIDX2*XX B 26 

C COMPUTE THE DIFFERENCE BETWEEN THE DATA AND THE B.G.R. B 27 

D1F=RH0(I )-BGR(I ) B 28 

VAL= VAL+D IF**2 B 29 

C TRANSFER AROUND THE UNWANTED GRADIENT CALCULATIONS WHEN JOKE IS 1 B 30 

IF ( JOKE • EQ. 1 ) GO TO 4 8 31 

D IF =2.#0 1 F 8 32 

C CALCULATE THE COMPONENTS OF THE GRADIENT (GRAD(l) ANO GR ADI 2 1 ) B 33 

GRAD ( 1 ) =GRAD (1)-DIF*TIDX2*XX B 34 

E XP 1 =E XP I X2D T 1 1 B 35 

E XP S= (E XPl+i . /E XP1— 2 . ) * T ( I ) 8 36 

GRA0(2)=GRAD(2) + {AK*BGRU > / X ( 2 > -X (1) / EXPS > *D I F 8 38 

4 CONTINUE B 39 

VLL= VAL B 40 

IF (J0KE.EQ.1I RETURN B 41 

DO 5 11=1,2 B 42 

5 GRDD (II) =GRAD 1 1 1 I 8 43 

RETURN 3 44 

END 3 45 - 


FCNCTION FUNKY (X) C 1 

C THIS FUNCTION CALCULATES THE INTEGRAND IN THE B.G.R. C 2 

COMMON /EXPK/ K C 3 

IF (X.EQ.O.) GO TO 1 C 4 

E XP 1 =E XP ( X) C 5 

FUNK Y=E XP1 +-1 • /EXP1— 2. C 6 

F UNK Y = X** K/F UNKY C 7 

RETURN C 8 

1 FUNKY =0. C 9 

RETURN C 10 

end c 11- 
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BLOCK DATA £ 1 

COMMON /ESTIM/ E ST ,E PS , LI MI T , I ER E 2 

COMMON /E XP K / K E 3 

C THE OROER OF THE B LOCH-GRUE NET SEN RELATIONSHIP MAY BE CHANGED. 

C 1U DO SO, CHANGE THE VALUE OF K IN THE FOLLOWING DATA STATEMENT. 

DATA Kit)/ E 4 

DATA EST »EPS » LIMI T/l. E— 2 , i . E— 5 ,1000/ E 5 

E ND E 6~ 


C 

C 

C 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

€ 

c 

c 

c 

c 

c 

c 

c 

c 


c 

€ 

c 

c 

c 

c 

c 

c 

c 

c 


12 


F 2 

SIRROITINE FLTPwL F 3 

F 4 

PURPOSE F 5 

TO FIND A LOCAL MINIMUM OF A FUNCTION UF SEVERAL VARIABLES F 6 

BY THE METHOD OF FLETCHER AND POWELL F 7 

F 8 

USAGE F 9 

CALL FLTPwLIFUNCT »N»X»F,G) F 1C 

F 11 

DESCRIPTION OF PARAMETERS F 12 

FUNCT - USER-WRITTEN SUBROUTINE CONCERNING THE FUNCTION TO F 13 

BE MINIMIZED. IT MUST BE OF THE FORM F 14 

SUBROUTINE FUNCT ( N ,ARG ,VAL , GRAD ) F 15 

AND MUST SERVE THE FOLLOWING PURPOSE F 16 

FUR EACH N-DI MENS I QNAL ARGUMENT VECTOR AR G, F 17 

FUNCTION VALUE AND GRADIENT VECTOR MUST BE COMPUTED F 18 

AND, ON RETURN, STORED IN V AL AND GRAD RESPECTIVELY F 19 

N - NUMBER OF VARIABLES F 20 

X - VECTOR OF DIMENSION N CONTAINING THE INITIAL F 21 

ARGUMENT WHERE THE ITERATION STARTS. ON RETURN, F 22 

X HOLDS THE ARGUMENT CORRESPONDING TO THE F 23 

COMPUTED MINIMUM FUNCTION VALUE F 24 

F - SINGLE VARIABLE CONTAINING THE MINIMUM FUNCTIUN F 25 

VALUE ON RETURN, I.E. F=F(X). F 26 

G - VECTOR OF DIMENSION N CONTAINING THE GRADIENT F 27 

VECTOR CCRRESPONDI NG TO THE MINIMUM ON RETURN, F 28 

I.E. G=G ( X ) . F 29 

fc ST - IS AN ESTIMATE OF THE MINIMUM FUNCTION VALUE. F 30 

EPS - TESTVALUE REPRESENTING THE EXPECTED ABSOLUTE ERROR. F 31 

A Rt ASON ABLE CHOICE IS 10**(-6), I.E. F 32 

SOMEWHAT GREATER THAN 10**(-D), WHERE D IS THE F 33 

NUMBER OF SIGNIFICANT DIGITS IN FLOATING POINT F 34 

REPRE SENTATI CN. F 35 

LIMIT - MAXIMUM NUMBER OF ITERATIONS. F 36 

IE R - ERRUR PARAMETER F 37 

IER = U MEANS CONVERGENCE WAS OBTAINED F 38 

IER = 1 MEANS NU CONVERGENCE IN LIMIT ITERATIONS F 39 

IER =-l MEANS ERRORS IN GRADIENT CALCULATION F 40 

IER = 2 MEANS LINEAR SEARCH TECHNIQUE INDICATES F 41 

IT IS LIKELY THAT THERE EXISTS NO MINIMUM. F 42 

H - WORKING STORAGE OF DIMENSION N*( N+7 )/2 . F 43 

F 44 

REMARKS F 45 

I) THE SUBROUTINE NAME REPLACING THE DUMMY ARGUMENT FUNCT F 46 

MUST BE DECLARED AS EXTERNAL IN THE CALLING PROGRAM. F 47 



C II) IEK IS SET TC 2 IF , STEPPING IN ONE OF THE COMPUTED F 48 

C DIRECTIONS, THE FUNCTION WILL NEVER INCREASE WITHIN F 49 

C A TOLERABLE RANGE OF ARGUMENT. F 50 

C IER = 2 MAY OCCUR ALSO IF THE INTERVAL WHERE F F 51 

C INCREASES IS SMALL AND THE INITIAL ARGJMENT WAS F 52 

C RELATIVELY FAR AwAY FRCM THE MINIMUM SJCH THAT THE F 53 

C MINIMUM WAS OVERLEAPED. THIS IS DUE TO THE SEARCH F 54 

C TECHNIQUE WHICH DCU6LES THE STEPSIZE UNTIL A POINT F 55 

C IS FOUND WHERE THE FUNCTION INCREASES. F 56 

C F 57 

C SUBROUTINES AND FUNCTION SUBPROGRAMS REQU IRED F 58 

C FUNCT F 59 

C F 60 

C METHOD F 61 

C THE METHOD IS OESCRIBED IN THE FCLLOWING ARTICLE F 62 

C R. FLETCHER AND M.J.O. POwELL, A RAPID DESCENT METHOD FOR F 63 

C MINIMIZATION, f 64 

C COMPUTER JOURNAL VCL.6 , ISS. 2, 1963, PP. 163-168. F 65 


C THIS SUBROUTINE IS A MODIFICATION OF THE FMFP PROGRAM FROM THE IBM F 66 

C SCIENTIFIC SUBROUTINE PACKAGE f a? 


C 


c 

c 

c 

c 

c 

c 


1 


2 

3 

4 
C 
C 

5 
C 

c 


c 

c 


SUBROUTINE FLTPWL (F UNCT » N ,X ,F ,G) 

COMMON /ESTIM/ EST, EPS, LIMIT, IER 
COMMON /SPACE/ HID 

DIMENSIONED DUMMY VARIABLES 
DIMENSION X ( 1 ) , G (1 ) 

COMPUTE FUNCTION VALUE AND GRADIENT VECTOR FOR INITIAL ARGUMENT 
CALL FUNCT IN,X,F ,G) 

RESET ITERATION COUNTER AND GENERATE IDENTITY MATRIX 
IER=0 
KOUN T =0 
N2=N+N 
N3=N2+N 
N 31 =113 + 1 

K=N31 

DO 4 J =1 , N 
H(K ) =1. 

N J=N- J 

IF (NJ) 5,5,2 
DO 3 L =1 , N J 
K L=K + L 
H(KL) =U. 

K =KL + 1 

START ITERATION LOOP 
KUUNT=KUUNT+1 

SAVE FUNCTION VALUE, ARGUMENT VECTOR AND GRADIENT VECTOR 
OLDF =F 
DO 9 J =1 , N 
K =N + J 
HIK) =G< J) 

K=K + N 
H (K ) =X ( J ) 

DETERMINE DIRECTION VECTOR H 
K=J+N3 
T =0. 


P 70 
IF 7i 
F 72 
F 73 
F 74 
F 75 
F 76 
F 11 
F 78 
F 79 
f 80 
F 81 
F 82 
F 83 
f 84 
F 85 
F 8 6 
F 8? 
F 88 
F 89 
F 90 
F 91 
F 92 
F 93 
F 94 
F 95 
F 96 
F 97 
F 98 
F 99 
F 100 
F 101 
F 1C2 
F 1C3 
F 1C4 
F 1C 5 
f 106 
F 107 
F 1C8 
F 109 

f no 
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6 


? 

8 

9 

c 


c 

c 

c 


10 

c 

c 

c 

c 

c 

0 

II 

c 

c 

€ 

c 


c 

c 

c 

13 

14 

15 
C 
C 

16 

C 

C 

17 

C 

c 

c 

c 

la 


DO 3 L-l.N 
I =I-G t LKHIK) 

IF { L” J ) 6 i 1 , 7 

K=K+N-L 
GOTO 8 
K ~K+ 1 
CON FI NUE 
H ( J J = T 

CHECK WHETHER FUNCTION WILL DECREASE STEPPING ALONG H. 

D Y ■" 0 © 

H NR M =0o 
G NR M “ G » 

CALCULATE DIRECTIONAL DERIVATIVE AND TESTVALUES FOR DIRECTION 
VtCTOR H AND GRADIENT VECTOR G. 

DO 10 J = i,N 

HNRM =HNRM+A6 S (H ( J } ) 

GNRM =GNRM+AB S ( G ( J) ) 

0 Y=DY+H ( JJ*G < J> 

REPEAT SEARCH IN DIRECTION OF STEEPEST DESCENT IF DIRECTIONAL 
DERIVATIVE APPEARS TC BE POSITIVE OR ZERO. 

If (DVI 11,54,54 

REPEAT SEARCH IN DIRECTION OF STEEPEST CESCENT IF DIRECTION 
VECTOR H IS SMALL COMPARED TO GRAOIENT VECTOR G. 

IF ( HNRM/GNRM-EPS) 54,54,12 

SEARCH MINIMUM ALONG DIRECTION H 

SEARCH ALONG H FOR POSITIVE DIRECTIONAL DERIVATIVE 

F Y=F 

ALFA=2.*(EST-F) /DY 

AMBD A =1 • 

USE ESTIMATE FOR STEPSIZE ONLY IF IT IS POSITIVE AND LESS THAN 
1. OTHERWISE TAKE 1. AS STEPSIZE 
IF iALFAS 15,15,13 
IF ( ALFA— AMBD AI 14,15,15 
AMBO A =ALF A 
A L F A "* G ® 

SAVE FUNCTION AND DERIVATIVE VALUES FOR OLD ARGUMENT 

F X=P Y 
D X=D Y 


STEP ARGUMENT ALONG H 
DO 17 I =1 »N 
XII ) =Xt I )+AMBDA*H( l ) 

COMPUTE FUNCTION VALUE AND GRADIENT FOR NEW ARGUMENT 
CALL FCNCT (N»X,F,G) 

F Y=F 

COMPUTE DIRECTIONAL DERIVATIVE DY FOR NEW ARGUMENT. TERMINATE 
SEARCH, IF DY IS POSITIVE. IF DY IS ZERO THE MINIMUM IS FOUND 

D Y — C . 

DU 18 I =1 ,N 
D Y=0 Y + GII i*H (I) 

IF ( 0 V » 19,39,22 


F 111 
F 112 
F 113 
F 114 
F 115 
F 116 
F 117 
F 118 
F 119 
F 120 
F 121 
F 122 
F 123 
F 124 
F 125 
F 126 
F 127 
F 128 
F 129 
F 130 
F 131 
F 132 
F 133 
F 134 
F 135 
F 136 
F 137 
F 138 
F 139 
F 140 
F 141 
F 142 
F 143 
F 144 
F 145 
F 146 
F 147 
F 148 
F 149 
F 150 
F 151 
F 152 
F 153 
F 154 
F 155 
F 156 
F 157 
F 158 
F 159 
F 160 
F 161 
F 162 
F 163 
F 164 
F 165 
F 166 
F 167 
F 168 
F 169 
F 170 
F 171 
F 172 


C 



c 

TERMINATE SEARCH ALSC IF THE FUNCTION VALUE INDICATES THAT 

F 

173 

c 

A MINIMUM HAS 8EEN PASSED 

F 

174 

19 

IF (FV-FXI 2U ,22,22 

F 

175 

c 


F 

176 

c 

REPEAT SEARCH AND DOUBLE STEPSIZE FOR FURTHER SEARCHES 

F 

177 

21 

AMBD A =AMBDA<- A LF A 

F 

178 


ALFA =AM6DA 

F 

179 

c 

END OF SEARCH LOCP 

F 

1 80 

c 


F 

181 

c 

TERMINATE IF THE CHANGE IN ARGUMENT GETS VERY LARGE 

F 

182 


IF <HNRM*AMBDA-1. ElO ) 16,16,21 

F 

183 

c 


F 

1 84 

c 

LINEAR SEARCH TECHNICUE INDICATES THAT NO MINIMUM EXISTS 

F 

185 

cl 

IER = 2 

F 

186 


RETURN 

F 

18? 

c 


F 

188 

c 

INTERPOLATE CUBICALLY IN THE INTERVAL DEFINED BY THE SEARCH 

F 

189 

c 

ABOVE ANU COMPUTE THE ARGUMENT X FOR WHICH THE INTERPOLATION 

f 

190 

c 

POLYNOMIAL IS MINIMIZED 

F 

191 

22 

T=0. 

F 

192 

23 

IF (AMBDAJ 34,39,24 

f 

193 

2** 

Z-3.*IFX-FY) /AMBDA+-OX+U Y 

F 

194 


ALFA=AMAX1 (ABS(Z) , AB S ( D X ) , ABS ( DY ) ) 

F 

195 


OALFA=Z/ALFA 

F 

196 


DALE A=DALF A*DALF A-DX/ALFA*D Y/ALFA 

F 

197 


IF (DALFA) 54,25,25 

F 

198 

2b 

W=ALFA*SQRT(DALFAi 

F 

199 


ALF A =D Y-0 X+ W+ W 

F 

200 


IF (ALFA) 26,27,26 

F 

201 

26 

ALF A = ( D Y- Z+ W ) /ALFA 

F 

202 


GO TO 28 

F 

203 

27 

ALFA = ( Z+D Y— w ) / 1 Z* DZ+- Z+ D Y) 

F 

204 

28 

ALFA=ALFA*AMBDA 

F 

205 


DO 29 1=1, N 

F 

206 

29 

XII )=X(I H-IT-ALFAJ *H (I) 

F 

207 

C 


F 

208 

C 

TERMINATE, IF THE VALUE OF THE ACTUAL FUNCTION AT X IS LESS 

F 

2C9 

C 

THAN THE FUNCTION VALUES AT THE INTERVAL ENDS . OTHERWISE REDUCE 

F 

210 

c 

THE INTERVAL BY CHOCSING ONE END-POINT EQUAL TO X AND REPEAT 

F 

2 1 i 

c 

THE INTERPOLATION. WHICH END-POINT IS CHOOS EN DEPENDS ON THE 

F 

212 

c 

VALUE OF THE FUNCTION AND ITS GRADIENT AT X 

F 

213 

c 


F 

214 


CALL FLNCT (N,X,F,G) 

F 

215 


IF (F-FX) 30,30,31 

F 

216 

36 

IF (F-FY) 39,39,31 

F 

217 

31 

DALFA =0. 

F 

218 


DO 32 1=1, N 

F 

219 

32 

DALFA=DALFA + G (I) *H ( I ) 

F 

220 


IF (DALFA) 33,36,36 

F 

221 

33 

IF (F-FX) 35,34,36 

f 

222 

34 

IF (DX-DALFA) 35 ,39,35 

F 

223 

35 

F X=F 

F 

224 


D X=D ALF A 

F 

225 


T =ALFA 

F 

226 


AM8DA =ALF A 

F 

22? 


GO TU 23 

F 

228 

36 

IF (FY-F) 33,37,38 

F 

229 

37 

IF (DY-OaLFA) 38,39,38 

F 

230 

38 

F Y=F 

F 

231 


D Y=DALFA 

F 

232 


AMBDA=AMBUA-ALFA 

F 

233 


GO TO 22 

F 

234 
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C F 235 

C TERMINATE, IF FUNCTION HAS NOT DECREASED DURING LAST ITERATION F 236 

39 IF (ULDF — F+EPS) 54,40,40 F237 

C F 238 

C COMPUTE DIFFERENCE VECTORS OF ARGUMENT AND GRADIENT FROM F 239 

C IRQ CONSECUTIVE ITERATIONS F 240 

40 Du 41 J = i ,N F 241 

K=N+J F 242 

l-UKl =G( Jl-HIK) F 24 3 

K=n+K F 244 

41 h C K ) = X i 3 ) — H IK) F 245 

C TEST LENGTH OF ARGUMENT DIFFERENCE VECTOR AND DIRECTION VECTOR F 246 

C IF AT LEAST N ITERATIONS HAVE BEEN EXECUTED., TERMINATE, IF F 247 

C BOTH ARE LESS THAN EPS F 248 

i £R =0 F 249 

IF ! KCUNT-N) 45,42,42 F 250 

42 T=u. F 251 

Z=0, F 252 

OQ 43 3=1 ,N F 253 

K=N+J F 2 54 

» =H s K I F 2 5 5 

K=K«-N F 256 

T=T+ABS(H (KJ ) F 2 57 

43 Z=Z+W*HIK) F 258 

IF ( HNRM-EPS) 44,44,45 F 259 

44 IF ( T-E P SI 59,59,45 F 260 

C F 261 

C TERMINATE, IF NUMBER OF ITERATIONS WOULD EXCEED LIMIT F 262 

45 IF < KOUNT-LI MI T) 46,53,53 F 263 

C F 264 

C PREPARE UPDATING OF MATRIX H F 265 

46 ALFA=0. F 266 

DU 5(j J = 1,N F267 

K=J+N3 F 268 

W=0» F 269 

DO 49 L = 1 ,N F 270 

KL=N+-L F 271 

W = rt+H < KL ) *H ( K ) F 272 

if ( L- J ) 47,48,48 F 273 

47 K=K+N-L F 274 

GO TO 49 F 275 

48 K =K + 1 F 276 

49 CONTINUE F 277 

K=N + J F 278 

A LF A =A LF A + I K ) F 279 

50 Ml J ) = W F 280 

C F 281 

C REPEAT SEARCH IN DIRECTION OF STEEPEST CESCENT IF RESULTS F 282 

C ARE NUT SATISFACTORY F 283 

IF { l* ALFA) 51 ,1,51 F 284 

C F 285 

C UPDATE MATRIX H F 286 

5 1 K =N J- 1 F 2 8 7 

DU 52 L =1 , N F288 

KL=N2+L F 289 

DO 52 3 = L ,N F 290 

NJ=N2’i-J F 291 

H(K) =H<K)+HI KL) *H (NJI/Z-H <L) *H< JI/ALFA F 292 

52 K=K-s-i F 293 

GO TO 5 F 294 

C END OF ITERATION LOOP F 295 

C F 296 
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c 

NO CONVERGENCE AFTER LIMIT ITERATIONS 

f 

29 7 

53 

I ER =i 

P 

298 


RETUkN 

F 

299 

C 


F 

300 

C 

RESTORE OLO VALUES OF FUNCTICN AND ARGUMENTS 

f 

301 

54 

DO 55 J=1,N 

F 

3C2 


K=N2+J 

F 

303 

55 

X ( J ) =H { K) 

F 

304 


CALL FUNCT (N,X,F ,G ) 

f 

3C5 

C 


F 

306 

C 

REPEAT SEARCH IN OIRECTICN OF STEEPEST DESCENT IF DERIVATIVE 

F 

30 1 

C 

FAILS TO BE SUFFICIENTLY SMALL 

F 

308 


IF (GNRM-EPS) 58,58,56 

r 

309 

C 


f 

310 

c 

TEST FOR REPEATED FAILURE OF ITERATION 

F 

311 

56 

IF ( I ER ) 59,57,57 

F 

312 

57 

I ER =- 1 

F 

313 


GO TO 1 

F 

314 

58 

I ER =0 

F 

315 

59 

RETURN 

F 

316 


END 

F 

317- 
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APPENDIX C 


EXAMPLES OF PROGRAM OUTPUT 
Example 1 - White and Woods Data of 1959 


SOURCE 

TEMPERATURE 

RHO 

RHO CALCULATED 

DIFFERENCE 

WHW059 

10.000000 

0.003200 

0.001020 

0.002180 

WHW059 

15.000000 

0.017000 

0.007719 

0.009281 

WHW359 

20.000000 

0.051000 

0.031382 

0.019618 

VJHW059 

25.000000 

0.120000 

0.036670 

0.033330 

WHW059 

30.000000 

0.230000 

0.133421 

0.046579 

WHW059 

AO. 000000 

0.540000 

0.499982 

0.040018 

WHW059 

50.000000 

0.950000 

0.936325 

0.013675 

WHW059 

60.000000 

1.430000 

1.436732 

-0.005732 

WHW059 

70.000000 

1.960000 

1.954969 

-0.004959 

WHW059 

80.000000 

2.500000 

2.501859 

-0.001859 

WHW059 

90.000000 

3.030000 

3.038230 

-0.008230 

WHW059 

100.000000 

3.550000 

3.570093 

-0.020093 

WHW059 

120.000000 

4.600000 

4.615721 

-0.015721 

WHW059 

140 .000000 

5.600000 

5.637964 

-0.037964 

UHW059 

160.000000 

6.650000 

6.640883 

0.009117 

WHW059 

180.000000 

7.650000 

7.628653 

0.021347 

WHW0S9 

200.000000 

8.600000 

8.504656 

-0.004655 

WHW059 

220.000000 

9.600000 

9.571464 

0.028536 

WHW059 

250.000000 

11.000000 

11.008536 

-0.008535 

A 

THETA 

SUM OF (0IFFERENCES)**2 



39.950176 

217.537354 

0.009226 





Example 2 - Cox Data of 

1943 



SOURCE 

TEMPERATURE 

RHO 

RHO CALCULATED 

DIFFERENCE 

MC0X43 

77.33000 0 

2.460000 


2.454481 

-0.004481 

MC0X43 

273.200001 

12.410000 


12.389841 

0.020159 

MC0X43 

373.400002 

17.180000 


17.193873 

-0.013873 

A 

THETA 

SUM OF (DIFFERENCES) **2 



39.512384 

210.770609 

0.000619 
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Example 3 - Cox Data of 1943 and White and Woods Data of 1959 


SOURCE 

TEMPERATURE 

RHO 

RHO CALCULATED 

DIFFERENCE 

WHW059 

10.000000 

0.003200 

0.000929 

0.002271 

WHW059 

15.000000 

0.017000 

0.007037 

0,009963 

WHW059 

20.000000 

0.051000 

0.028769 

0.022231 

WHW059 

25.000000 

0.120000 

0.030298 

0.039702 

WHW059 

30.000000 

0.230000 

0.171999 

0.358001 

WHW059 

40.000000 

0.540000 

0.478832 

0,051168 

WHW059 

50.000000 

0.950000 

0.910253 

0.039747 

WHW059 

60.000000 

1.430000 

1.411158 

0,318342 

WHW059 

70.000000 

1.960000 

1.943838 

0.016152 

WHW059 

80.000000 

2.500000 

2.437599 

0,012401 

WHW059 

90.000000 

3.030000 

3.032195 

-0.002195 

WHW059 

100.000000 

3.550000 

3.572969 

-0.022959 

WHW059 

120.000000 

4.600000 

4.637077 

-0. 037077 

WHW059 

140.000000 

5.600000 

5.677648 

-0,077643 

WHW059 

160.000000 

6.650000 

6.698330 

-0.048330 

WHW059 

180.000000 

7.650000 

7.703241 

-0.353241 

WHW059 

200.000000 

8.600000 

8.695818 

- 0 . 0 9 5 8 1 B 

WHW059 

220.000000 

9.600000 

9.578715 

-0.078715 

WHW059 

250.000000 

11.000000 

11.139177 

-0,139177 

MC0X43 

77.330000 

2.450000 

2.342040 

0, 117953 

MC0X43 

273.200001 

12.410000 

12.259842 

0.153158 

MC0X43 

373.400002 

17.180000 

17.045341 

0.134659 


A THETA SUM OF ( DI F FERENCE S ) **2 

41.621598 223.468170 0.114080 
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TABLE I. - COX DATA OF 1943 


Temperature, 

Total electrical 

K 

resistivity, 


P» 


^ohm-cm 

77. 33 

2.46 

273.2 

12.41 

373.4 

17. 18 


TABLE II. - WHITE AND WOODS 


TABLE III. - BLOCH-GRUENEISEN CONSTANTS AS FUNCTION 


DATA OF 1959 


OF MAXIMUM TEMPERATURE 


Temperature, 

K 

Total electrical 
resistivity , 

P, 

pohm-cm 


Maximum temperature, 
K 

Constant of metal, 

A, 

pohm-cm 

Characteristic temperature 
of resistivity, 

e R’ 

K 

10 

0.0032 


50 

27.730 

1.85 . 664 

15 

.017 


60 

31.233 

194. 136 

20 

.051 


70 

34.525 

202.274 

25 

. 12 


80 

36.568 

207.447 

30 

.23 


90 

37.530 

209.995 

40 

.54 


100 

37.905 

211.035 

50 

.95 


120 

38.588 

213.063 

60 

1.43 


140 

38.726 

213.490 

70 

1.96 


160 

39.372 

215.580 

80 

2. 50 


180 

39.777 

216.932 

90 

3.03 


200 

39.813 

217.056 

100 

3. 55 


220 

40.010 

217. 758 

120 

4.60 


250 

39.950 

217. 537 

140 

5. 60 




160 

6. 65 




180 

7. 65 




200 

8.60 




220 

9. 6 




250 

11.0 




273 

12. 1 




295 

13. 1 
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Figure 1. - Fletcher-Powell algorithm. 













Figure 2. - Main program. 
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Total electrical resistivity, p, pohm-cm 





ature of resistivity, 





Temperature, K 


Figure 6. - Constant of metal as function of highest temperature used. 



Temperature, K 

Figure 7. - Characteristic temperature of resistivity as function of maximum temperature used. 
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